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ABSTRACT 


Reexamining a suggestion of Cox and Smith (1974), it is found that 

intersecting supernova remnants can indeed generate and maintain hot 

-2 -3 6 

interstellar regions with n ^ 10 cm and T ™ 10 K- These regions 

are likely to occupy at least 30% of the volume of a spiral arm near 

the midplane of the gaseous disk if the local supernova rate there is 
, -7 -1 -3 

greater than 1.5 x 10 Myr pc . Their presence in the interstellar 
medium is supported by observations of the soft X-ray background. 

The theory required to build a numerical simulation of inter- 
acting supernova remnants is developed. The hot cavities within a pop- 
ulation of remnants will become connected, with varying ease and speed, 
for a variety of assumed conditions in the outer shells of old remnants 
Rayleigh-Taylor instabilities are important in forming connections. 
Apparently neither radiative cooling nor thermal conduction in a lai‘ge- 
scale galactic magnetic field can destroy hot cavity regions, if they 
grow, faster than they are reheated by supernova shock waves, but inter 
stellar mass motions disrupt the contiguity of extensive cavities 
necessary for the dispersal of these shocks over a wide volume. 

Monte Carlo simulations show that a qua si- equilibrium is reached 

7 

in the test space within 10 yrs of the first supernova and is charac- 
terised by an average cavity filling fraction of the interstellar 
volume. Aspects of this equilibrium are discussed for a range of super 
nova rates. Two predictions of Cox and Smith are not confirmed within 
this range: critical growth of hot regions to encompass the entire 
medium, and the efficient quenching of a remnant's expansion by inter- 
action with other cavities. 
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HOT INTERSTELLAR TUNNELS: I. SIMULATION OF 


INTERACTING SUPERNOVA REMNANTS 
I. INTRODUCTION 

It lias been suggested by Cox and Smith (1974; hereafter CS) that intersec- 
tions of supernova remnants contribute to the production of gas at about 10 6 K 
which is apparently observed in the interstellar medium. Evidence for this gas 
was first found at low galactic latitudes in soft X-ray emission exceeding the 
extrapolated diffuse background determined at higher energies (e.g. Bunner, 
et al. 1972). Kraushaar (1973) and Williamson, et al. (1974) have concluded 
that the soft X-ray observations require widespread thermal emission regions 
interwoven with cooler gas, because of the short mean free path for 250 eV 
photons in interstellar gas (Brown and Gould 1970) and because a broad range 
of other emission mechanisms can be eliminated. The analysis of ultrasoft 
X-ray observations by Levine, et al. (1976) is also consistent with this hypothesis. 
Since Vanderhill, et al. (1975) and others have set fairly stringent upper limits 
on the space densities of soft X-ray stars, the diffuse soft background seems to 
require hot interstellar gas. Contrary to early expectations (CS), ultraviolet 
observations of O VI interstellar absorption lines in the spectra of unreddened 
hot stars, e.g. Jenkins and Meloy (1974) and York (1974), give no unambiguous 
evidence for this gas, since the observed absorption may be associated with the 
stars themselves (Castor, McCray, and Weaver 1975). However, the latter 
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work does not weaken the X-ray evidence, since hot gas produced by these uncom- 
monly hot stars (earlier than ~B2; Snow and Morton 1976) cannot supply the 
diffuse soft X-ray background (McCray 1976). Moreover, Scott (1975) has 
recently given substance to the suggestion in CS that very hot regions may be 
preferred channels for cosmic ray propagation, and he points out that fresh 
estimates of cosmic ray lifetimes in the Galaxy favor this idea, although he 
oversimplified the problem. 

These considerations justify a continued exploration of the ideas of CS. 

The fundamental prediction of the CS model is that intersections of supernova 
remnants (SNRs), which involve the interaction of older low-density cavities 
left by SNRs with strong shock fronts of other supernovae (SNe), generate more 
hot gas from an initially cold interstellar medium (ISM) than in the absence of 
interactions. A high SN rate per unit volume or a population of very-large- 
diameter SNRs might generate enough hot gas to satisfy the X-ray observational 
demand by simple superposition, but it turns out that overlaps, interactions, 
and cooperative regeneration of remnants are still likely to occur. In this paper 
the important physical processes encountered in the model are discussed. 

A fairly detailed numerical simulation of interacting SNRs has been devel- 
oped, and a description of its design and results occupy the latter half of the 
paper. The complex geometry and competing processes which connect and 
rejuvenate remnants can only be treated by numerical simulation except in the 
most trivial cases. We have investigated the large-scale three dimensional 
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behavior of a test section of the gaseous galactic disk under the influence of 
evolving and interacting SNRs. The accuracy of the method rests on the validity 
of models for individual and collective SNR behavior. Ultimately the models 
must come from separate gas-dynamical calculations, although even when such 
results are available the complexity of the simulation forces simplification of 
the models. Thus the simulation should provide reliable global properties of 
the medium but not detailed information about small-scale local structure. The 
chief aim of the discussion here will be to show that SNR intersections can 
quickly generate large volumes or "tunnels” of very hot gas from a cold start- 
ing medium, under conservative and reasonable assumptions. The two most 
important ones, explained later, are that SNRs can be treated as spheres and 
that mass motions may be treated implicitly. 

In the next few sections the SNR models selected as building blocks for the 
simulation are assembled. First, in § II models for a noninteracting, "isolated" 
SNR and the ambient medium are discussed. §m examines pairs of interacting 
remnants and the mechanism by which their central cavities can become con- 
nected, while §IV deals with the evolution of larger aggregates of remnants. 
After these developments, the simulation will be described in §V and its 
results presented and discussed in §VI. Conclusions are summarized in the 
final section. 
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ISOLATED SUPERNOVA REMNANTS 


a) Very Late Remnant Evolution 

Tlie SNR radius of interest for the SNR connection mechanism of GS is not 
the shock radius R s but the inner radius R ; which delimits the hot, low-density 
central cavity (see Cox 1972 for an overview of SNR evolution). In old remnants 
Rj is the inner radius of the dense shell which forms around the cavity. The 
cavity must be reached by another SNR shock before a connection can proceed. 
There are many uncertainties in the later evolution of R. , but the particular 
evolution chosen here is sufficient to portray the general behavior of all solu- 
tions with respect to connections. 

Spherically symmetric models such as Model A of Chevalier (1974) can be 
extrapolated in time with some confidence to the point where the interior pres- 
sure P t (which is well defined because of the high sound speed in the cavity) 
equals the ambient pressure P ISM , under the assumption that the dense shell 
remains hydrodynamically stable. This point is normally reached by an age of 
10 6 yrs, when the final radius R^ is in the range 35 - 70 pc, depending weakly 
on the initial blast energy E Q ( < 10 5 1 ergs) and the ambient density n Q 
( > 0.1 cm' 3 ). Due to the residual momentum of the dense shell, the cavity 

may actually expand beyond the pressure equilibrium point. Then P. would 
drop below 3? ISM , since the rate at which edge material can be accelerated 
inward by external (interstellar) pressure is inadequate to enforce pressure 
equilibrium. In the present work isolated remnant sizes will be restricted to 
the smaller and simpler regime. 
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The representative SNR evolution which will be used exclusively for the 
remainder of this paper is given in Table 1. It is Chevalier’s (1974) Model A 
scaled up to E 0 = 6x10 5 0 ergs and n Q « 0.7 cm' 3 by the method described in 
his paper. For each value of the remnant age t 5 = t/10 5 yrs, the shock radius 
and speed and the fraction of E 0 remaining as thermal energy in the hot cavity, 

Ej /E 0 , are given in columns 2 through 4. Miese quantities were taken from 
Chevalier’s Figure 1 as far as possible and then extrapolated using Ms power 
law fits and assuming the later expansion of the cavity to be adiabatic. Assum- 
ing only that the mean pressure in the dense shell (which forms just before 
t = 1.0) is proportional to the ram pressure of material entering the shock, 
one can calculate the mean density compression factor in the shell <x> = 
< n sheii> / n o 1 c °l umn 5. Since essentially all the SNR mass then resides 
in the shell, mass conservation gives the shell thickness AR and implies an 
inner radius R. =R S - AR. Finally the cavity pressure is P. - E. /BttR. 3 . 

Talcing account of the continued accretion of cavity gas onto the inside of 
the dense shell (which does not cause significant departures from the adiabatic 
drop of Ej ), one finds a cavity density n. < 10' 2 n Q at R. f . With this density 
at P gM , the isobaric cooling time is quite long. It turns out that a faster 
"destruction” of the contiguous, low-density region within a SNR may occur 
when it is pinched or bridged and becomes disconnected into separate volumes 
by motions of cooler gas. Therefore a cavity may "die” for purposes of the CS 
model after a time of roughly 
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Myr. 


(1) 


T 


SNR 



To simplify interpretation of the simulation results, r SN R - 4 Myr has always 
been used. But if SNR connections become frequent, some cavities are never 
destroyed, in the sense that the low-density regions never completely disappear. 


b) Cloud Structure 

Realistic models of the ISM should take its density inhomogeneities into 
consideration in a self-consistent way. For example, the evolution of SNRs 
depends qualitatively on the cloud structure of the ambient medium. Clouds 
which are small in the dimension parallel to a shock velocity, as well as all 
sufficiently low-mass clouds, are accelerated by the SNR shock passage and the 
postshock flow, becoming part of the kinetic energy field of the remnant and 
eventually joining the dense shell. In contrast, massive or monolithic clouds 
left behind the shock do not seriously affect its dynamics (McKee and Cowie 
1975, Sgro 1975). Thin filamentary or lamellate clouds are probably responsi- 
ble for some of the 21-cm emission and absorption (Heiles 1973, Verschuur 
1974; cf. Greisen 1976, Spitzer and Jenkins 1976). In a study of conditions in 
hot regions Shapiro and Field (1976) derive a high pressure which implies that 
some H I clouds used to obtain P 1SM could be lamellate rather than spheroidal. 
In fact, tlnn cloudlets arise inevitably from abundant interstellar shocks due to 
SNe themselves, as well as stellar winds, nova outbursts, etc. Elmergreen 
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(1975, 1976) argues that the ionized intercloud density is of order 1 cm' 3 . A 
consistent pietuie of the cycle of gas, "evaporated" from clouds by stellar 
ultraviolet to form (clumpy) intercloud material and then swept up by shocks 
into new clouds, would include thin clouds as a prominent stage because of 
thermal instabilities behind shocks (Chevalier and Theys 1975; McCray, Stein, 
and Kafatos 1975). Therefore if one is practically constrained to uniform- 
density SNR models (such as Chevalier 1974), one should use a value of n 
somewhere between the average and intercloud densities, as we have here. 

Note that McKee and Cowie claim that only the filling fraction of clouds is 
important for evaluating their dynamical effect on SNRs: We contend that the 
density which can be quickly accelerated by the shock is the significant quantity, 
and that it depends on cloud morphology. 

A significant consequence of the flattening of clouds swept up by shock 
fronts is that the embedded magnetic fields tend to be compressed and aligned 
chiefly parallel to the filaments and sheets. Heiles (1973) has discussed obser- 
vations which support this expectation. This alignment affects transport prop- 
erties between clouds and hot regions, and it may also mean that a magnetic 
field contrast can be maintained. 

c) Thermal Conduction 

The conduction of heat from hot, low-density regions to surrounding or 
included cooler regions probably does not critically affect the survival of the 
hot regions. Let us consider two possible effects in turn. 
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Does a hot, old SNR cavity cool by losing heat to its dense shell? The 
thermal conductivity in the absence of a magnetic field for T ~ 106 K is 
(Spits er 1962) 

k 0 ~ 6 X 10* 7 T it2 ergs' 1 cm' 1 K' 1 . (2) 


The conductivity perpendicular to a magnetic field is decreased by 


3 ~ 10 


l B 


r-3 


(3) 


where n p is the proton density and B is in /jlG. If the dense shell of an old 
remnant has expanded in the presence of a large-scale field component B Q 
which is uniform over scales larger than a remnant, then the shell will have a 
parallel field component everywhere except at two small "caps." When conduc- 
tion to the shell is predominantly perpendicular to the shell field, solution of the 
heat flow equation (e.g. Spitz er 1956) with shows that heat transport in this 
situation is utterly negligible. The unit area of the shell where frozen-in field 
lines might directly connect the hot interior and the cold shell gas can be esti- 
mated from flux-freezing: 


cap 

TiR . 2 


B. 


rn. 


2/3 


~ — ~ [ — ] < 10' 4/3 


B, 


(4) 
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P 'er these small areas, heat flow is one dimensional and along the field, so the 
heat flow equation can be solved for the energy flux, using Spitzer's (1956) 
method. Mass flow is neglected, pressure equilibrium at 10~ 12 dyn cm* 2 is 
assumed, and a recent radiative cooling coefficient* is used. For a typical 
cavity temperature of 10 6 K the quasi-steady conduction region is ~ 11 pc 
thick and the conduction cooling time is 


t 


cond 


19 

N 



(5) 


where E. is in ergs and A„„„ in pc 2 . N is the number of unit exchange areas 
over which field-free conduction is established; N - 2 for a spherical SNR with 
"caps" but no included clouds threaded on the cavity field. 

Does "evaporation" of clouds destroy cavities? Chevalier (1975) claims 
that evaporation of "standard" clouds will cool a cavity in the total absence of 
magnetic fields in < 10 6 yrs. Since fields are probably present, (4) estimates 
the upper limit to a unit field-free exchange area for a standard cloud, unless 
better estimates become available from detailed studies of the heat flow prob- 
lem such as Cowxe and McKee (1976). These authors have reached conclusions 
different from Chevalier’s, finding that saturation of conductivity and the sensi- 


*L(T) 1.2 X 10* 19 T-" 2 eig cm 2 s"* is an excellent approximation for 3 X 10 5 <T <5 X 10 2 K to the total cooling curve 

of Raymond, Cox, and .Imitli (1976) and recovers the same temperature dependence as that deduced by Spitzer (1956). 


tivity of cloud evaporation rates to the hot gas temperature prevent even 
embedded "standard" clouds unshielded by magnetic fields from cooling a hot 
region much below 10 6 K. Moreover, thin clouds may be much less affected by 
conduction. In view of this controversy, we may appeal to soft X-ray observa- 
tions as suggested by Chevalier (1975) to reassure us that considerable gas at 
~ 10 6 K does survive conduction cooling. However, note also that the arguments 
of McKee and O striker (1975) suggest that only by admitting some cloud evap- 
oration can one keep interacting SNR cavities from thermal runaway toward 
higher temperatures ; new SNRs occurring in low-density regions have higher 
cavity temperatures than the first SNRs in a cold medium, 

HI. CONNECTIONS BETWEEN CAVITIES 

The probability P^ that a new SN occurring at random will be inside an old 
remnant cavity is equal to the fraction f of the volume occupied by cavities at 
that time. If SNRs do not intersect, f is equal to 1 - e' q , where q is the dimen- 
sionless porosity defined in CS: 

q “ rr SNR v snr » ( 6 ) 

where r is the SN rate per unit volume and V gN R = 4 7rR jf 3 /3. For large q, 
intersections are inevitable, and f ^ 1 - e' q . Porosity has a simple meaning 
for q « 1, so that intersections are rare, namely f = p in « q. Those SNe 
which do not occur inside cavities may yet occur near enough to a cavity to 
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burst into it. Let us consider how cavities are rejuvenated and/or connected 
for these two types of interactions, i.e, for SNe inside and outside cavities. 

The outcome of a SN event inside a remnant cavity is a shock expanding 
rapidly through low-density gas out to the edge; cf. Chevalier's (1974) model 
with n 0 = 10 ' 2 cm" 3 or Cowie (1976). The old cavity will be strongly reheated, 
and the resulting general repines surization would give a renewed impetus to an 
enclosing dense shell. If cold inclusions are present in the cavity, the strong 
shock will go around them to reach all the low-density volume, and it will drive 
into them as into the dense rim of the cavity. The rejuvenation of the cavity is 
due to reheating of cavity gas and to compression of dense inclusions. Clouds 
struck by shocks have been considered by McKee and Cowie (1975) and Sgro 
(1975), McKee and Cowie concentrate on the "main" cloud shock structure. 
Sgro's numerical results show that the shock enters a dense cloud from all 
sides, and tends to keep it compressed if the postshock cloud gas radiatively 
cools before the cloud is completely shocked. Less dense clouds reexpand 
sooner. The critical cloud density is given, using Sgro's notation but the more 
recent L(T) given in §H, by nj /n ~ 10* 3 (nd) - 1 ^ 2/3 v s 4 / 3 , where n is the 

low cavity density and n c ' is the corresponding critical value of the cloud 
density. The size of the cloud is d pc, v s is in Ion s _I , and 0 > 1 is the 
ratio of postshock pressures inside and outside the cloud, discussed both by 
Sgro and by McKee and Cowie. For n ~ 10* 2 , d ~ 5, and v J s = 1000, inclusions 
withn c > 1cm* 3 should remain relatively compressed after the shock passes. 
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With respect to the cavity demise represented by equation (X), shocks forestall 
the filling, bridging, and pinching off of cavities by clouds; since the act .al life- 
time thereby added to a cavity is difficult to pin down quantitively, a range of 
rejuvenation hypotheses should be studied. 

If a SN occurs outside a cavity but nearby, and if the strong shock front of 
the young SNR considerably overlaps the cavity, at least one low-density channel 
may be opened between the two cavities. Such channels, forming a "connection,” 
ensure that a later shock propagating through one cavity will propagate through 
the connection into the other, providing a rejuvenating effect similar to that 
outlined above. The establishment of the connection takes some time, during 
which the newer SNR expands normally, but if its interior pressure is diverted 
through a channel before its shell forms, this remnant may not reach its full 
normal size. 

How do these channels come about? Dense gas in an old SNR shell which 
is struck from the outside by a shock is first accelerated suddenly as the front 
passes, and then further accelerated and decelerated in a complex way which 
depends on the density structure through the shell. When the shock speeds up 
into lower density on the far side of the shell, a rarefaction wave propagates 
back through the shell. With pressure differences across the shell of more 
than factors of 10 2 on top of the various accelerations, it is plausible that the 
shell is subject to Rayleigh-Taylor instabilities. Full calculations, which take 
into account the effect of compressibility emphasized by Chevalier and Theys 
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(1975), are quite difficult to perform for these instabilities, but we can estimate 
the magnitude of the connection time delay using the approximate treatment of 
Frieman (1954). This assumes a static, incompressible, nonviscous fluid sub- 
ject to an acceleration. Idealize the situation as follows: When the new shock 
hits the old shell, the shock speed drops suddenly and a reflected shock forms, 
increasing the thermal pressure between the shocks. If the average shock speed 
in the shell is v , then the transmitted shock takes a time t ~ AR/v to 

cross the shell. For the important case in which the younger active SNR can 
still be described by a Sedov adiabatic blast wave and has no shell, Table 1 gives 
^ cross ~ 10 5 F rs * 'The thermal pressure difference A? across the shell is 

sufficient to cause an instability after t cross has elapsed and P, of the new 
remnant has meanwhile dropped from overall expansion. This AP implies an 
effective gravity g ~ AP / p AR, where p is the average shell density and 
p s AR is constant during the compression of the shell by the crossing shock. 

A rough lower bound for g under conditions favorable to connection is 4 km s" 1 
(10 5 yr)' 1 . The dispersion relation for the instability can be written (Frieman 
1954) X~2tt gt RT 2 , where t R T is the growth time in the linear regime and 
X is the scale (both parallel and perpendicular to the shell) over which the 
instability has grown aftpr an elapsed time t . As an indication that the shell 
ruptures we require X~AR, calculating a growth time ~ 10 5 yrs. Com- 
paring t cross + t R T to time scales for SNR evolution, and recognizing that 
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the two delays overlap somewhat if the instability sets in during t cross , one 
can conservatively conclude that som s , SNR intersections result in connections. 

However, when a SNR is struck from the outside, it may not be necessary 
to break through a uniform shell as thick as that suggested by AR in Table 1. 
The above instability develops from small perturbations in the shell, but one in 
fact expects large perturbations to be present in SNR shells. When a dense 
shell forms behind the spherical shock front of an isolated SNR, the shock speed 
drops suddenly and the cooling gas following the shock is strongly decelerated 
as it runs into the slowing shell. It is very likely that a Rayleigh-Taylor insta- 
bility grows at this stage (Kahn 1975). If the preshock medium is cloudy, the 
instability takes on the character of the Chevalier-Theys mechanism. Unfortu- 
nately the magnetic field topology has an important influence on the detailed 
nature of the instability. The dynamical effects of tangled and uniform fields in 
the shell are compared briefly by Chevalier (1974). Nonuniforraities in isolated 
SNR shells also arise after about 10 6 yrs in a uniform interstellar field 
because flux tubes bent by the SNR expansion, tend to straighten and move 
toward the explosion center on a time scale given roughly by the radius of 
curvature of a bent flux tube divided by the Alfven speed in the tube. 

Clearly the corrugated, filamentary, or fragmentary shells resulting from 
the preceding processes will develop low-density connection channels more 
quickly than uniform shells when struck by strong shocks. SNR cavity connec- 
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tions seem virtually assured, although a proof by hydrodynamical calculation 
has only been attempted by Jones (1975a), whose computer code was inadequate 
to provide a conclusive test. 

IV. EVOLUTION OF MULTIPLE CAVITIES 

A well established ch: 'n or cluster of cavities is similar to an isolated 
cavity in having hot, low-density gas which is cooling slowly and low magnetic 
field strength. The boundaries or "walls' 1 of the cavity are composed of portions 
of old shells of component SNRs and, in part, of clouds or "inclusions" which 
were enveloped by SNR shock fronts or have migrated into the cavity. The hot 
cavity gas reacts to at. least three processes: 1) Distortion and spatial redis- 
tribution by interstellar motions (CS), 2) cooling and thermal conduction, and 
3) sporadic reheating and repressurization above P ISM by fast shock waves 
reverb exiting through the cavity. Consider the third process. 

a) Shock Waves in Tunnels 

When a strongly expanding SNR intersects a chain (if it is very long we may 
call it a "tunnel"), it may eventually initiate a strong shock wave in the low-density 
gas, promptly for inside SNe and after some delay for outside SNe. In one likely 
topology, tunnels may be considered as having low-density branches leading 
away from the locality of the intersection, so that the shock divides into "tunnel 
shocks” into various branches and "wall shocks" into high-density matter. The 


wall shocks involve more and more wall surface as the faster tunnel shocks 
advance. 

The shock fronts generated in a tunnel by one supernova, along with their 
enclosed volume, will be called a shock system. Figure 1 shows a simple shock 
system involving a cylindrical tunnel. The rate of change of the volume of this 
shock system is 


TV 

dt 


= 2 A + S v, 


St 


(7) 


where A is the cross-sectional area of a branch of the tunnel undisturbed by the 
new shock system, S is the wall surface area already affected by the system, 
and the factor of 2 refers to the two branches in Figure 1. The tunnel shock 
speed v st and the wall expansion speed v are derived in the Appendix, and 
their different character is related to the different sorts of pressure jumps at 
tunnel shocks and wall system boundaries. Tunnel shocks are approximated by 
pressure waves with a jump P 2 - , where P refers to the inshocked 

tunnel and P 2 to the shock system. When P t is comparable to P, , it does not 
cost much energy to rejuvenate the tunnel gas, and the tunnel shocks may be 
weak although they are still fast. One may assume that tunnel shocks add volume 
and perhaps some energy (from P ) to the system without doing much work. 
However, the wall expansion constitutes an essentially adiabatic movement 
against P IS M which does work, the work appearing as power radiated behind the 
slow shock. Therefore 
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dP 2 
dt 


\_ 

V 


(P 2 -P 1 )2Av st + j P 2 Sv v 


( 8 ) 


Dividing (8) by (7) yields an equation for dP 2 /dV which has two simple limiting 
solutions: In the early phase of a shock system v st is large, P., »P ^ and A 
is still comparable to S. Then one can neglect terras like Sv w and find P 2 oc V* 1 , 
But after a time, Sv w becomes larger than Av s[ , and volume is being added to 
the shock system by wall expansion much more rapidly than by tunnel shock 
propagation. Then P 2 V' 5 ^ 3 until it approaches P ISM and levels off because 
v w c s , the sound speed in the walls. In a tunnel of arbitrary geometry the 


terms containing A should be replaced by a sum over all the branches of the 
tunnel, but a time is usually reached for a strong shock system when wall expan- 
sion dominates dP 2 / dV. 


b) An Example 

It is appropriate now to present an example of SNB, connection and of a 
two-remnant or binary shock system. For convenience the example will be 
discussed within the framework of the computer simulation, which is not 
described until the next section, and the opportunity is taken to suggest the 
modest level of accuracy with which a simulation can represent physical processes. 
Consider the representative evolution of a pair of SNBs shown in Figure 2. The 
second one occurs nine 50,000-year time steps after the first and 60 pc away. 

The pressures and radii would roughly follow Table 1 until both the interaction 
and the connection delay t crQss + t have occurred. In Figure 2 this delay 
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begins at the "inter section” time t Jnt and lasts until the "connection” time 
t con * After t , P. must drop on the average more rapidly than in the isolated 
remnant evolution, due to simultaneous expansion of the shock system into the 
old cavity and spherically. To ensure that the time -aver aged system pressure 
is not overestimated in the simulation where the actual hydrodynamics are not 
done explicitly, the two pressure drops are calculated separately. On each 
time step after t con , a trial integration is performed of the pressure drop AP cav 
which would occur if the shock had started across the old cavity at t con driven 
by the current P 2 * P; of the shock system. The current P, is lower than at 
t , because the new SNR continues to evolve normally. The trial propagation 
beginning at t con has a duration At p (which in the particular case of Figure 
2 is larger than a time step). If the result cf the trial indicates that the propa- 
gation would be at least half over within the current time step, then the cavity 
shock takes effect as if it occurred instantaneously at this time, t cav . Other- 
wise the only drop in P, for that step is due to normal spherical expansion (if 
a binary system, or adiabatic wall expansion over the current surface area S 
if the system already covers more than one SNR). At the chosen t cav , the pres- 
sure drop AP cav is well approximated by the result of the calculation shown by the 
dashed line in Figure 2a, AP cav is subtracted from P 2 , which thereafter 
describes the pressure of both cavities, since pressure equilibrium is roughly 
maintained throughout a shock system (see Appendix). The compression of 
cavity shock propagation into a single time step is required by the simulation, 
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since there is no way to represent the position of a shock which has crossed a 
cavity partially; see §V. The final linear increase in radii in Figure 2b is due 
to the relaxation of P 2 to P JS by expansion over the wall surface S as mentioned 
above. 

Note that the most effective shock systems for rejuvenating tunnels are 
actually initiated by SNe within a cavity rather than outside as in the example 
above. In such stronger shock systems, of which Figure 1 is a contrived exam- 
ple, there is no connection delay and the first cavity shock is completed when P 2 
is still quite high, perhaps enabling several cavities to be shocked in quick suc- 
cession. 

c) Lifetimes of Tunnels 

The rejuvenating influence of shocks on the cavities they pass through has 
been demonstrated in § §11 and IE. Here it may be added that for two cavities 
connected by a previous binary interaction, a later tunnel shock tends to reopen 
the low-density channels that form the essential connection: A shock compresses 
inclusions which clog the channels, and the divergence of the shocks leaving the 
connection tends to move the walls away from the channels if the Chevalier- 
Theys mechanism is operative. 

A small amount of reheating of the cavity gas every few million years keeps 
ahead of radiative cooling. Thermal conduction in a very old section of tunnel 
could become important if the heat exchange area perpendicular to the field 


increases substantially with time. Otherwise equation (5) tells us that the 
replenishment of 2 x 10 4 9 ergs every 19/N million years years to each unit 
cavity volume { ~ IQ 6 1 cm 3 ) keeps ahead of conduction. If the effect of both 
cooling and conduction is to decrease the cavity pressure below P [SM , a shock 
repressurization to, say, 10% above P ISM both pays the energy debt and fore- 
stalls the collapse of the cavity. Such a repressurization is quite easily supplied 
by tunnel shocks in the models to be discussed. Most tunnel shocks provide at 
least 1.3 P, SM . 

Of course if rejuvenating shocks fail to appear for r SNR , equation (1), the 
channels and tunnels may become clogged. If a cavity is allowed to cool, it is 
probably first filled by intercloud gas at 10 3-4 K from evaporation of clouds by 
the ultraviolet flux of hot stars which will be free in cavities; cf. Elmergreen 
(1976). But many regions survive through rejuvenation so long that they are 
removed from the disk principally by buoyancy as discussed by Jones (1973). 
This process will not be considered in detail in this paper, but it should be 
remembered that the buoyant force on a hot "bubble” is zero at the midplane 
and increases with height as long as the Z-eomponent of the galactic gravity 
increases. Therefore at the midplane where the SN rate per unit volume is 
highest and tunnels are most prevalent, buoyancy is weak. 
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V. METHOD OF SIMULATION 


Model SNe are distributed randomly in space and time throughout a test 
volume representing a square section through the gaseous galactic disk in a 
spiral arm in the solar neighborhood. As time passes one must calculate and 
maintain information about the lifetimes, intersections, and shock systems of 
the collection of active and inactive living SNRs. Special care is needed to make 
conservative estimates of imperfectly understood phenomena and the somewhat 
arbitrary parameters that arise in the equations; where errors are likely, we 
have tried to err on the side that minimizes tunnel growth. 

If a numerical experiment is to be meaningful, it is necessary to specify 
carefully which of the possible experiments has been performed. Several 
aspects of the simulation should be mentioned at the outset. First, cavities are 
spherical, being represented by a birth location and time and by a radius R . 
However, a calculational method has been found to allow for tangled or cloudy 
tunnel structure when propagating tunnel shocks. Secondly, no mass motions of 
the gas are followed explicitly. Motions are implicitly treated in the models 
already described for evolution of isolated remnants and propagation of shocks 
in chains of cavities. Thirdly, SNRs can occur only on the mesh points (spacing 
10 pc) and time steps (50,000 yrs). The resulting granularity is small compared 
with significant distances and time scales for SNRs and is on the order of uncer- 
tainties which are generated by subsequent manipulations of SNRs present in the 
test space. 
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The simplest version of the simulation which has been investigated assumes 
that SNRs expand instantaneously to a radius of 40 pc, chat any two overlapping 
living SNRs instantly make a perfect connection, and that one new intersection 
with a chain of remnants instantly rejuvenates the entire chain, no matter how 
long. With these unreasonable assumptions, the growth of a tunnel system 
depends solely on the value of the porosity q, equation (6). However, for porosi- 
ties so low that tunnels did not grow, they cannot grow for any assumptions. 

The results indicate q must be greater than about 0.01, whereas galactic SNRs 
certainly amount to more than this. The second version of the simulation, the 
results of which were quoted in CS, contained a considerably more realistic 
treatment, but there were errors and unnecessary assumptions which have been 
corrected in the version described in this paper. Smith (1975) describes these 
differences in detail, but from now on only the latest version will be considered. 

The size of the space represented was X X Y X Z = 400 X 400 X 250 pc 3 . In 
order to prevent spurious surface effects at the artificial X and Y boundaries, re- 
peating boundary conditions were imposed to mimic the behavior of the larger 
galactic disk; this means for example that an expanding SNR having a center at low 
enough X (or Y) can expand across the boundary and make a connection with a 
SNR at high X{or Y). In contrast, surface effects on the top and bottom of the 
test space will occur and must be allowed, because the thickness of the real disk 
is not much larger than the SNR. diameter. Remnants near these boundaries 
are allowed to expand into the volume above the test space, representing a 

galactic halo, although no SNe occur outside the space. 
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The Z -distributions of SNe and gas have not been treated comprehensively 
at this stage of the work. That of SNe has been included, but the gas distribu- 
tion is modeled by the disk with constant density n Q sandwiched between two halo 
atmospheres with a lower constant density n H , High -altitude SNRs will inter- 
act with the hot, flattened halo which various authors have proposed and which 
is predicted on the basis of the CS model. Such a halo is supplied with hot 
diffuse gas by convection from the disk, e.g. Jones (1973) and with mechanical 
energy by shock waves of high -altitude SNRs, e.g. Chevalier and Gardner (1974). 
The SN distribution in Z was taken to be a decreasing exponential with scale 
height 90 pc, symmetric about the center of the test space (Z = 0) as proposed 
by Ilovaisky and Lequeux (1372). This distribution is applied here only out to 
Z = 120 pc where it is effectively truncated. Higher-altitude SNRs will affect 
only the halo, whereas we are primarily interested in the disk. 

In the time steps following the birth of each SNR, it was expanded according 
to Table 1. If a SN occurred within an existing cavity, its entire blast energy 
was available to begin a shock system. This approximation is not too bad since 
the time scale for radiation is longer for an explosion in a lower density (Che- 
valier 1974). Lower fractions of E Q are available for systems resulting from 
intersections made only after some expansion. For increased resolution during 
a remnant T s initial step when almost half the expansion takes place, tills step 
was divided in half as indicated in Table 1. To remain conservative in the matter 
of connections, an overlap of the spheres of radii R. of two SNRs was required 
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to initiate a connection, rather than an overlap of one shock radius R g and one 
inner cavity radius R. . If the calculated delay t + T found at t. was 
so long for an intersecting pair that negligible pressure would have been avail- 
able after the delay had elapsed, the connection was rejected and the two SNRs 
remained unlinked. In the borderline case, a viable connection was made but 
the older SNR was not rejuvenated; a later shock system could propagate through 
such a connection. 

Shock systems arising from SNe in cavities and from successful connections 
were expanded by the method of §IV. In doing this, it was necessary that the 
entire volume AV of a SNR which is added to the system be shocked in one 
uninterrupted procedure; otherwise the shock fronts would normally have been 
left delineating partial SNR volumes, and some further specification of the 
spatial extent of SNRs would have been necessary beyond their centers and radii. 
The maintenance of identifiable component SNRs in a chain is a weakness of the 
method but is also what makes it feasible. In complex situations where a 
shock system simultaneously expands into branches with fast local tunnel shock 
speeds and branches with slow speeds, the order in time of cavity shocks and 
pressure drops was carefully observed. 

The tunnel shock area A in equations (7) and (8) was evaluated in the geom- 
etry of intersecting spheres. The value of S to be used is a more important 
quantity, since wall shocks predominate for most of the duration of a system. 
One might continue to visualize a tunnel as a set of overlapping spheres, but 
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this would seriously underestimate the surfaee-to-volume ratio S/V which a 
region presents to a propagating shock system. In the computer code, old 
tunnels appear tangled on dimensions smaller than isolated SNR sizes, i.e. 
similar to embedded cloud sizes , in spite of the restricted spatial information, 
by letting each SNR have an S/V which increases with time t after the 
SN blast time t gN according to the law 


t - 1 


i?(t) = i? f + U, - V exp 


SN 


(9) 


'sv 


Here the initial ratio is that of a cylinder of radius 30 pc, i.e, 77 . = 0.067 pc 1 , 
while the final ratio attained is that of a cylinder of radius 5 pc, i.e. 
i? f “ 0.4 pc * 1 . The age of the SNR, t - t , may be many isolated remnant 
lifetimes. The time constant r sv has been given the value 4 x 10 6 yrs to 
agree with the lifetime hypothesis. Almost any tangling assumption has an 
enormous effect on the distance a shock system can extend in an old section of 
tunnel. To obtain AS swept into the shock system, as a volume AV is swept 
out, one calculates AS = r\ (t) AV, using r\ for the SNR into which the 
shock propagates. 

SNe occurring within 30 pc of the test space boundaries in the Z-direction 
(the top and bottom three layers of mesh points) were assumed to make "con- 
nections" with the diffuse halo whether or not connections were made with other 
SNRs. In keeping with the assumption of uniform density, a halo “connected 
SNR which became part of a shock system was assumed to be the center of a 
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hemispherical halo shock front which began after the shock system had rejuve- 
nated the SNR and continued with increasing radius into the halo, contributing 
to the drain on shock system pressure until P 2 -*■ P TSM . A more accurate 
treatment might consider the actual gas distribution at high Z and the lion- 
spherical nature of high-altitude SNRs (Chevalier and Gardner 1974). The halo 
shock speed was given by equation (A.5) using P (I = 8 x 10' 13 dyn cm" 2 
instead of P and p H - 6 x 10‘ 27 g cm' 3 instead of ; these values are 
suggested by models of a halo arising from hot material welling up from the 
tunnels (Cox 1974; cf, Spitzer 1956) or from accretion of local intergalactic 
matter (e.g. Cox and Smith 1976). This very approximate view of the disk and 
halo can be compared with theoretical models of the Z- distribution of matter 
such as Kellman’s (1972). However, the aim here was simply to include some 

i 

representation of the existing sink for high-altitude SN energy and avoid its 
spurious concentration to the disk. 

If a shock system managed to shock all the SNRs in a short chain or cluster 
before depressurization to P ISM , it was then allowed to continue pushing out its 
walls until equilibrium with P 1SM was reached. This adiabatic expansion affects 
only the outer walls and not the inclusions. Tlius the surface S used in inte- 
grating (7) and (8) for this purpose is just the outer boundary of some configur- 
ation of overlapping spheres. The expansion produces a AR t during one time 
step which is added to each Rj in the system, and can increase R { above 
R. f =* 46.7 pc. An arbitrary upper limit of 60 pc was placed on R ; , because 
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without buoyancy one remnant could expand to fill the test volume, an unrealis- 
tic situation which would be an artifact of the assumption that SNRs retain their 
identity within a chain. 

The rejuvenation of a SNR by a shock system was accomplished by revising 
its projected death time from t SN + r SNR (or its most recent value) to equal that 
of the SNR which initiated the rejuvenating shock system. However, its life 
expectancy was never reduced thereby. This estimate of the degree of rejuve- 
nation reflects the fact that an older shock system cannot extend a SNR's life- 
time as much as a newer system, but as mentioned previously the exact amount 
of the extension must be regarded as unknown. Since the degree of rejuvenation 
is a critical consideration, one comparison run with a high SN rate was made 
in which a remnant's death time could be extended upon rejuvenation by only 
half the above amount. As discussed in the next section, the fraction of the 
whole test volume ending up in tunnels was somewhat lowered by this change. 

At the end of each time step, all isolated SNR that had lived to their limit 
were totally erased from the test space, and their positions effectively restored 
to uniform n Q wherever living remnants did not extend. All chains that had 
one or more members due to die were analyzed to erase only those SNRs which 
were too old; members do not have the same death time unless they were most 
recently rejuvenated by the same shock system. For example, if a long, straight 
chain were very old when it was intersected by new SNRs at both ends, the 


27 


middle of the chain might die and be erased before the two shock systems met, 
leaving two pieces no longer having a connection. 

Considerably more detailed simulations can certainly be imagined, but for 
generalized chains even this approximate formulation is already rather 
complicated to encode for the computer, 

VI. RESULTS AND DISCUSSIONS 

Results will now be presented for a sequence of simulations having the 
same values of E 0 , n Q , and r SNR , but different SN rates per unit volume, r. 
Table 2 lists the parameters used in the various runs. The fraction <f > 
of the test volume contained In SNR cavities (angle brackets denote averages 
over Z) is the most basic calculated quantity and as used in the following 
includes both isolated and connected SNRs. To give the reader a feeling for the 
fluctuations in <f > for different r, <f> is plotted in Figure 3 as a function 
of time for the runs, along with the "random input function." The same random 
set (taken from Rand 1955) has been used for all runs, merely varying the rate 
at which SNe are put into the space, so that in each run they appear at the same 
positions and in the same order and relative spacing in time. Although the runs 
are short, the poor statistics are compensated by the clarity with which the 
effect of varying r can be isolated. To test for bias of the common data set, 
which would affect all runs, two additional data sets were compared for one 
rate. No bias tending to increase <f > was found. 
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Although the results depend separately on <r> , t snr , and V SNRJ one may 
characterize the runs by their values of <q> since only <r> is varied. 

The exponential Z “distribution of SNe has been truncated in such a way that up 
to 26% of the: total galactic SNe would fall outside the space, depending on how 
far the distribution holds. Thus in Table 2 the Z “averaged r is given for each 
model, and the equivalent galactic rates are derived for truncation at Z - 120 pc 
and Z = . 

Let us now focus individually on the types of results obtainable from this 
study. 

a) Death and Rejuvenation of Remnants 

Table 3 presents the distributions of SNR ages in the form of percentages 
of SNRs with ages in the ranges given in the first column. Only SNRs which 
died by the end of the run are included. The average age at death increases 
steadily with <q> for the sequence of Models A, B, and C. In fact, a short 
run performed for <q> - 0.30, which is not shown in Figure 3, gave an aver- 
age age more than twice as long as r SNR . These average values appear on 
Figure 3 as " r cff 

b) Shock Waves in Tunnels 

Tables 4a and b give the distributions of peak cavity shock speed v (e.g. 
at the highest point on the dashed line in Figure 2a) and of the peak wall speed 
v w , for each unit cavity which was crossed by a shock system throughout the 
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run. During At prop both speeds drop. The v sl distribution is peaked at c g 
but has a second maximum at the high end due to strong shocks from SNe within 
cavities, and the v w histograms have a similar structure. 

In Table 4c are propagation times At prop for the same sample of shocks. 

The breadth of the distributions, which cover an order of magnitude, is due to 
differences in At prop between fast v st ’s across small cavity volumes AV and 
slowv^’s across large AV's. 

c) Fraction of Interstellar Volume within Tunnels 

The duration in time represented by the runs In Table 2 and Figure 3 is 
about 15 <q> _1 Myrs, Within this time each run attained a quasi-steady 
state in which <fl> averaged more than <q> . If SNRs enjoyed no cooper- 
ative interactions, one would obtain <f> = 1 - e" <c t > . Here it is apparent 
that the extension of SNR lifetimes through rejuvenation has added more cavity 
volume than has been lost from the incomplete expansion of SNRs which break 
into 'ther cavities. 

The quantity <f > averaged over time, < f > , is plotted in Figure 4 
against <q> with approximate uncertainties in the quasi-steady value esti- 
mated from Figure 3. Models A, 3, and C are barely consistent with a straight 
line through the origin, although a better fit is given by a function which is 
concave upward. Figure 4 includes the run at <q> = 0.30 mentioned previously. 

The critical value of q, which CS Implied would separate slow tunnel growth 
from rapid and extensive growth, was not found here for <q> < 0.30, 
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i.e. SN rates < (25 yr gal)* 1 . The tendency of all runs to reach a quasi- 
steady <f > does not encourage a search at higher q for such a critical value. 
However, by performing very long runs one might detect a critical q-dependence 
for a secular increase in <f> above these quasi -steady levels. 

A general property of tunnel growth in these runs, which summarizes the 
results so far and explains the lack of critical behavior, will now be discussed: 
Tunnels are primarily built and enlarged by intersections between very old 
SNRs and new ones in the Sedov phase, but once tunnels develop they are kept 
rejuvenated by SNe occurring within the.lr volume. Isolated SNRs spend most 
of their time fully expanded, and a matrix of binary interactions like the 
example in §IVb shows that such old SNRs with thick shells (R. f Z 47 pc, 

A R ~ 13 pc) are hard to rejuvenate. It can be done by a SN inside the cavity 
or by one outside but still in the Sedov phase, R ( < 24 pc, when it intersects. 
The probability that a SN will occur within a distance R of a very old SNR 
within its lifetime is p out = 1 - e' aq , where a = (R^ + R) 3 /R, f 3 . In CS we 
claimed that R = R jf would produce a connection, i.e. a= 8. From the re- 
quirement that the active SNR be in the Sedov phase we find a - 3.5. For <q> 
= .25, the mean probability for being rejuvenated once by an outside SN is thus 
.63, but the weak shock system usually stops after the single SNR. On the 
other hand, with <q> = .25 there is a fairly high SN rate within the tunnels 
themselves, each SN capable of rejuvenating V — 2 * 10 62 cm 3 (CS), i.e. to a 
distance R t ~ 200 pc, if such a fully connected tunnel network has been built 
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that V t is available. Once rejuvenated at some time, the volume is guaranteed 
a second rejuvenation by another occurrence within a time r SNR - t (where 
t prop ~ 2 Atprop for the chain). Since the average tunnel shock speed for Model 
C was about 350 km s’ 1 , one estimates t ~R /v ~ 0.6 Myr. Hence the 
probability of rejuvenating V t is p^ ~ 1 - exp [-r(T SNR - t prop )V { J ~ .97, 
in an extensive tunnel network. But SNBs originating within tunnels add little 
or no volume to the tunnel. The comparison of p out and p jn emphasizes the 
connection-damping effect of SNR shells, 

d) Z -Distribution of Volume Fraction 

Figure 5 shows a typical relationship, at the end of Model B, between the 
distribution of SNe with height above the plane, represented here by q(Z), and 
the resultant fraction f(Z) of the volume left in tunnels and SNRs. The values 
of f(Z) for positive and negative displacements from Z - 0 have been averaged 
together; the breadth of the low-Z peak is a measure of the fluctuation in the 
position of the maximum f(Z) around Z = 0, roughly ±20 pc. At high Z we see 
that f(Z) ~ q(Z) in agreement with the low effective porosity there and the low 
probability for SNR interaction. The curves for f(Z) and q(Z) separate near 
Z = 100 pc, and f(Z) substantially exceeds q(Z) near the midplane. The SN 
rate per unit volume there is — 1,7 times the Z -averaged rate, which causes 
tunnels to grow first in this vicinity. Over a very long time they might spread 
to higher Z, which would provide a secular increase in the average <f > . 

For all four runs, the relation of f(Q) to q(0) is almost the same 
as that between <f> and <q> in Figure 4. 
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The f(Z) vs. Z curve should become flatter for a model in which both the 
scale height of ambient gas and the density decrease with Z are taken into 
account, because SNRs with the same E 0 will reach larger R. f in lower n Q . In 
fact it is possible that a correct treatment of high-Z Stilts with the same Z- 
distribution used here would show f(Z) increasing with height, a question which 
deserves further study: The slope of f(Z) vs. Z depends critically on the com- 
petition between the scale heights of the gas and SN distributions, as suggested 
by Reynolds (1975). The interaction between tunnels and the galactic halo 
ultimately depends on similar considerations. The eventual outcome of the 
development of high-altitude tunnels is the expulsion of hot gas from the galactic 
disk (cf. Paper I, Spitzer 1956, Jones 1973, Shapiro and Field 1976). 

e) Radii of SNR Cavities 

In these simulation models the radius of a SNR depends on its history in 
spite of the fact that all begin identically. A very young remnant that drives a 
shock into a neighboring cavity is rapidly depressurized and may not reach its 
full radius; after t. nt in such a case the remnant expands until t cav when the 
first cavity is actually "shocked."* Furthermore, a SNR whose expansion has 
ceased may renew it if subsequently repressurized, unless Rj - GO pc. 

With these aspects of the simulation in mind we present the distributions 
of radii observed for the runs in Table 5. The average SNR radius is less than 

*Without in chiding tliis expansion after intersection it is possible to find parameters for which the simulation gives f < q 
due to strong quenching of SNR expansion. 
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roughly R if foi* all runs and decreases with increasing q. A similar decrease 
is observed for the average radius at intersection. A comparison of Tables 5b 
and c shows the substantial increases of Rj between t jm and t cay . 

Hie results for the radii when one takes a snapshot of the test space at a 
typical time may be summarized very roughly by constructing graphs of the 
number of SNBs with diameter D. less than a certain value. A set of such 
curves, each averaged over a run, is given in Figure 6. On the time scale of 
the simulation, millions of years, there are never an appreciable number of 
remnants with Dj < 50 pc. The simulation does not need to resolve such small 
diameters, since the evolution of SNRs at this stage is mostly unaffected by 
interactions. Most SNRs predicted by these runs will be larger and luminous in 
soft X-rays, not radio emission, because a cavity cooling time is required for 
the decay of the X-ray emission (cf. Naranan, et al. 1976), whereas the radio 
emission decays quickly into the background continuum. The soft X-ray back- 
ground excess is probably a SNR superposition of some sort, but SNRs cannot be 
distinguished individually nearly as far in X-rays as in radio emission because 
of absorption, so that we cannot begin to construct a catalog of large galactic 
SNRs. Nonetheless, it is interesting to compare in Figure 6 the curve of Clark 
and Caswell (1976) for radio SNRs. Several features of these curves viewed 
together are evident: First, the ’’radio" and "X-ray” curves have essentially 
no overlap; the turnover in the radio curve above D = 35 pc is probably a 
limiting-suzTace-brightness effect. Secondly, large-diameter SNRs in the runs. 


34 


Dj - 50 - 100 pc, have N(<D.) curves which are always steeper in general 
trend than the approximate D ; 4 proportionality expected for "momentum- con- 
serving" SNRs which have formed dense shells. This can be understood as 
incomplete expansion of interacting SNRs; many SNRs have smaller diameters 
on the average than would be expected if all expanded independently. This 
steepening is most evident for = 50 - 70 pc in Model C. Thirdly, if one 
scales the number of SNRs in the 400 x 400 pc 2 test space by the total area of 
the disk, accounting for a nommiform SN rate over the disk, one finds a plausible 
correspondence between the predicted radio and X-ray curves, e.g. as shown by 
the dashed line extrapolation in Figure 6. For this line, chosen arbitrarily to be 
determined by Model B at large a scale factor of 10 3 was assumed, which 
would make sense if one third of the galactic disk (R - 12 kpc) were subjected 
to a rate per unit volume equivalent to a uniform galactic rate of (40 - 50 yr)' ! . 
Maza and Van den Bergh (1976) have recently confirmed that SNe seem to be 
somewhat concentrated to spiral arms. Also, Ogelman and Maran (1975) have 
proposed that young stellar associations are sites of multiple, correlated SN 
explosions. Because of this possible clustering and because the growth of 
tunnels depends on the SN rate per unit volume, the above steepening may occur 
even with an average galactic rate as low as the (150 yr)* 1 favored by Clark and 
Caswell. There is a hint of steepening in their radio N{ <D) curve for D = 30 - 
35 pc which could possibly correspond to the steepening predicted by the tunnel 
model. It would then be predicted that the thousands of soft-X-ray-emitting 
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SNRs in the Galaxy should show a steep N( <Dj) curve in intermediate diameters, 
although it is not clear how many can be observed. Seward, et al. (1976) have 
compared radio and X-ray fluxes from SNRs and find that many are unexpectedly 
faint in X-rays. However the larger distance scale adopted by Clark and Caswell 
leads to many disagreeing^' with Seward, et al. on the distances and diameters 
of SNRs and tends to weaken tius conclusion. 

f) Effective Remnant Volumes 

The volume of a SNR that remains isolated is simply V SNR =* 4.3 x IQ 5 pc 3 , 
but if it initiates a shock system, then its effective volume is the eventual volume 
of that system. Both cases are included in Table 6 which lists V cff for SNRs 
which died and for shock systems which were erased by the end of each run. 

A broad range of values is found, particularly for the higher SN rates. Note 
the high values in the second bin indicating that a common class of shock systems 
propagate through only about one cavity or twice V SNR , in agreement with the 
discussion at the end of §VIc. 

VH. CONCLUDING REMARKS 

Dense, thick shells insulating SNR cavities provide effective damping on 
the rate of remnant connections, but they by no means prevent all SNR inter- 
actions. Remnants of 6 x 10 50 erg explosions can break through shells <: 10 pc 
thick if they begin while still in the Sedov phase, because a Rayleigh-Taylor 
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instability occurs which is enhanced by irregularities in the shell. SNRs which 
occur outside in relatively dense locations and break into a tunnel add volume 
to the tunnel network but initiate only weak shock systems, while SNe which 
occur inside tunnels add no volume (taking the shape of their surroundings) 
but create very strong shock systems. Lifetimes of cavities, both single and 
multiple, probably depend chiefly on magnetic field topology over scales smaller 
than a remnant, since thermal conduction cooling depends sensitively on field 
geometry and radiative cooling does not destroy cavities. 

Numerous assumptions must be made to make progress on a global simula- 
tion. The resulting tunnel filling fraction is strongly influenced by sizes of 
isolated SNRs, which depend in turn on the choice of E 0 , , and the point at 

which expansion is halted (P. -»■ P JSM ). Cloudy inclusions within tunnels can be 
accounted for roughly by using a growing surface-to-volume ratio for each SNR 
cavity. The combination of a uniform density disk, sandwiched between uniform 
lower-density haloes, with a realistic Z -distribution of (identical) SNRs gives 
tunnel filling fractions which fall off with Z. Within this conceptual framework 
the numbers of effectively free parameters is understandably high, and the 
parameter space has not been exhaustively explored. However, fixing one’s 
choices for most of them, the dependence on the SN rate per unit volume as 
independent variable can be understood. For the life and death of the popula- 
tion of galactic SNRs as a whole the chief dependences are on V gKR , r SNR , and 
the degree of rejuvenation resulting from shocks which propagate through a 
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cavity. Critical growth of the tunnel filling fraction f for certain minimum values 
of porosity q (CS) was not found here, but rather a quasi -equilibrium in which 
f q. The quenching of SNB expansion by early interactions was found to be 
less effective than claimed by CS and by Jones (1975b). 

The presence of a local tunnel network in some region, e.g. the solar 
neighborhood, depends on the SN rate per unit volume in that vicinity over the 
preceding IQ 7 yrs, which may easily exceed the galactic average. SNe are 
concentrated at low Z (e.g. Clark and Caswell 1976), near spiral arms (Maza and 
Van den Bergh 1976), and perhaps near each other (Ogelman and Maran 1976). 
With respect to spiral structure, the density wave theory predicts passage 
through an arm in — 10 8 yrs, so that cunnels maj' be quickly established there 
and then decay in the interarm region. If galactic shocks are suppressed in a 
10 6 K gas with high sound speed, various observations may require this interarm 
tunnel decay (Chevalier 1976). Soft X-ray spectral studies (Burstein 1976, 
Naranan, et ah 1976, Levine, et ah 1976) are still unable to prove a thermal 
origin for the emission* or to detect whether the solar system is near to or 
surrounded by hot SNR cavities which could account for all the observed soft 

I 

i 

I X-ray background excess. 

i 

i 

I We do not accept the narrow view of Clark and Caswell (1976) that all 

) 

: galactic SNRs are known, because SNRs should be X-ray emitters long after 


♦Since the emission from tunnels is overwhelmingly collisionally excited lines (Kato 1976, Raymond, Smith, and Cox 1976, 
Shapiro and Moore 1976), the resolution of soft X-ray lines would be decisive. 
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they cease to be bright radio objects. But if their radio results are a good 
indication of present SN properties such as E 0 , n 0 , and the SN rate, then the 
structure of very hot interstellar regions may differ from the visualization of 
the present paper, which is based on earlier studies. Clark and Caswell find a 
ratio E 0 /n 0 ~ 5 x 10 51 erg cm 3 and a mean interval between SNe of < 150 yrs, 
both numbers being larger than found by earlier workers {e. g. Ilovaisky and 
Lequeux 1972). The final radii of such SNRs would be very large, and a rela- 
tively small number of such SNRs may in fact constitute sufficient cavity volume 
to explain the soft X-ray background excess, as suggested by Naranan, et al, 
(1976). In this case the porosity may still be large and SNR interactions thus 
frequent, but the scale of cavities is much larger compared to the thickness 
of the gaseous disk than in the models discussed in this paper. Interactions 
between SNRs would then assume less importance while those between SNRs 
and the galactic halo might determine the structure of hot interstellar regions, 
making fewer and generally wider cavities, rather than narrower but extensive 
tunnels. 
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APPENDIX 

Shock system expansion into tunnels and walls is of different 
character, as one can discover from the expansion speeds. The 
equations of conservation of mass, momentum, and energy across a 
shock front are 


Pl v l = R3 v 2 


(A.l) 


P l +Pl v l 


a _ 


= Po + 


P2 v 2 


( A .2) 


5kT, +^mv, 2 + — = 5k1o + •— mv 0 s + — (a -3) 

121 8nn L ^22 

where p is the mass density and m is the mass of an average nucleus. 
For tunnel shocks one notes that the magnetic field is negligible 
in tunnels, so that the equations yield a compression factor 


x 


= P2_ 
Pl 


4P 2 +R 


and a shock speed 

,4 




v 


st 


= (- 


p Ji.p i /s 


Pi 


(A. 4) 


(A- 5) 


Immediately after a SN within a tunnel the program might have 
P 1 ~ ^ISM an( * % ~ ^30 resulting in v st > 1500 Ion s" 1 with 

= 2xlQ _se g cm" 3 . More typical speeds are five to ten times lower 
(§VI). As the expansion of the shock system proceeds, P 2 rapidly 
decreases and v s{; drops to the sound speed, c g 90 km s _1 . If P^ 
is higher than P^g^ due to recent passage of a different shock system, 
v st is somewhat increased. 

At a shock front in the walls, the magnetic field ~ x s * 3p,G 
(where x s is the residual shell compression) is compressed anew along 
with the gas. These shocks are much slower than tunnel shocks due 
to higher p 1 , and radiative cooling is important just behind the 
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shock. Since Eg “ x in equation (a- 3), the wall shock compression 


- ^ ^ l/S 

* SW 


(A. 6) 


where P w ^ = S| S /8 u + P cr is the sum of magnetic and cosmic ray 
pressures. The shock speed is 


V — Vi 

sw V 1 


=/3l [i + (!Bi) 1/2 ]j 1/2 


(A. 7) 


^Pw 

However, as shown in the inset of Figure 1 , this is not the speed 
of the inner edge of a wall of the shock system, which is 


v w v sw 


Cl - “) 

x 


= {*-*- [l + tf£i) ] 

Pw E 


1 /a 'l l /a 


v , 1 / 2 

[i - (£» l) ]. 


(A. 8) 


For a mass density in the walls of p - 3xlCT s4 g cm" 3 , and the above 
value of P 2 , v ws =- 120 Ion s" 1 . The value of p w reflects the typical 
value of x g , about 2. Values of v w have a lower limit of the Alfven 
speed v^ ^ 6 1cm s" 1 . 

Because the diffuse tunnel gas is hot and has a high sound speed, 
the pressure in different parts of the shock system tends to equalise. 
Just behind a shock front the pressure will be higher than average; 
however, if one assumes a power law decrease in density behind the 
front and a velocity pattern which maintains this density structure, 
it can be shown (Cox 1974) that the kinetic energy par unit mass of 
the motion associated with a tunnel shock is less than one sixth of 
the thermal energy per unit mass of the material near the shock, for 
any ratio of preshock to postshock pressure /P 2 . Both kinetic 
energy and the driving pressure P 2 cause the shock to advance. Since 
the total energy content of the shock system is shared between kinetic 
and thermal energy, we overestimate P 2 somewhat by assuming that all 
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energy fed into the system is thermal and that the interior of the 
system is in pressure equilibrium. A workable approximation is to ignore 
the kinetic energy when calculating the decay of P 2 , because the lower 
than true rate of system expansion computed without this kinetic energy 
is increased by the compensating effect of larger P 2 . 
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Table 1 


Calculated SNR Evolution, n 0 = . 7 cm' 3 , E 0 = 6 x 10 5l) ergs. 


l 5 

(yrj 

(pc) 

V s , 

(km s' 1 ) 


<X> 

AR 

(pc) 

R , 

(pc) 

P. 

(10‘ 12 dyn cm' 2 ) 

.25 

17.6 

276 

.72 

4 

0 

17.6 

430 

.5 

22.2 

182 

.63 

4 

0 

22.2 

190 

1.0 

28.2 

84.3 

.17 

20 

0.5 

27.7 

25 

1.5 

31.5 

63.7 

.088 

14 

0.8 

30.7 

10 

2.0 

34.5 

52.3 

.074 

10 

1.2 

33.3 

6.5 

2.5 

36.9 

44.8 

.063 

8 

1.6 

35.3 

4.7 

3.0 

39.1 

39.5 

.057 

6 

2.3 

36.8 

3.7 

4.0 

42.7 

32.4 

.048 

4.0 

3.9 

38.8 

2.7 

5.0 

45.8 

27.8 

.043 

3.5 

4.9 

40.9 

2.1 

6.0 

48.4 

24.5 

.040 

3.0 

6.0 

42.4 

1.7 

7.0 

50.8 

22.0 

.038 

2.7 

7.1 

43.7 

1.5 

8.0 

53.0 

20.0 

.036 

2.5 

8.3 

44.7 

1.3 

10.0 

56.8 

17.2 

.034 

2.1 

10.8 

46.0 

1.1 

12.0 

60.1 

15.2 

.033 

1.9 

13.3 

46.8* 

1.0 

♦E-inal cavity radius R. f . 











Table 2 


Parameters Used in the Runs 



A 

Model 

B 

c 

1. * 

<q>, volume-averaged porosity 

.15 

.20 

.25 

.25 

Volume-averaged SN rate within 120 pc 

of Z = 6 (yr' 1 gal' 1 ) 

. . . 66* 1 

50* 1 

4i* ! 

41* 1 

SN rate including high-Z SNe; rate 

for whole galaxy (yr* 1 gal 1 ) 

49' 1 

37' 1 

3 0' 1 

30* 1 

<r>, SN rate per unit volume at 

1Z1 = 5 pc (10-7 Myr> pc-3> 

. . . 1.5 

2.0 

2.4 

2.4 

Elapsed time (Myr) 

. . . 98 

69 

52 

59 

SN blast energy E 0 (ergs) 

• • • 

6 X I0 50 



Ambient ion density n Q (cm* 3 ) 

. . . 

.7 



Isolated SNK lifetime r SNR (yr) 

. . . 

4 X 10 6 



Final isolated SNR volume V SNR (pc 3 ) 

• • • 

4.3 X 10 s 



“Tangling” time constant r sv (yr) 

. . . 

4 X I0 6 



Tunnel ion density p (g cm' 3 ) 

• • • 

2 X 10* 26 



Ambient interstellar pressure P |SM (dyn cm' 2 ) 

. . . 

1.0 X I0* 12 



Halo ion density p„ (gem* 3 ) 

. . . 

6 X IQ* 27 



Halo pressure I’u (dyn cm' 2 ) 

. . . 

8 X 10* 13 



•Half the rejuvenation lifetime extension of other runs, sec §§111 and V. 
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Table 3 

Distributions of SNR Ages at Death 


Model 

T eff Bin 

(Myr) A BCD 

4 (%) 57 45 43 34 

4-6 (%) 13 13 20 29 

6-8 (%) 12 16 7 18 

8-10 (%) 8 8 6 8 

10-12 (%) 6 5 9 4 

12- 16 (%) 2 5 5 3 

16-20(3) 1 4 8 2 

>201%) <1 3 3 <1 

Average r, ff (Myr) 5.8 7.2 7.5 6.4 
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Table 4 


Distributions of Quantities for Cavity Shocks 


Speed Bin 
(km s' 1 ) 


A 

Model 

B 

C 

I) 

a) PEAK TUNNEL SHOCK SPEEDS v jt 





< 100 

(%) 

13 

12 

14 

14 

100 - 200 

(%) 

54 

54 

54 

56 

200 - 300 

(%) 

8 

10 

12 

9 

300 - 400 

(70 

“> 

4 

3 

4 

400 - b00 

(70 

3 

3 

3 

3 

600 - 1000 

(70 


i 

1 

1 

> 1000 

(70 

18 

15 

14 

13 

Average v st (kin s' 1 ) 

b) PEAK WALL EXPANSION SPEEDS v w 

401 

360 

346 

332 

< 15 

(70 

63 

64 

62 

66 

15-30 

(70 

14 

17 

20 

17 

30 - 45 

(70 

3 

3 

4 

3 

45 - 75 

(70 

1 

1 

1 

1 

75 - 90 

(70 

8 

n 

i 

5 

6 

90-120 

(70 

1 

<i 

1 

<1 

120 - 135 

(%) 

9 

8 

8 

7 

Average v w (km s' 1 ) 

c) SHOCK PROPAGATION TIMES At n 

y P ro P 

Time Bin 
(10 5 yr) 

31 

28 

27 

26 

< 1 

(%) 

18 

19 

27 

17 

1 - 2 

(70 

25 

28 

24 

27 

2-3 

(70 

17 

21 

21 

23 

3-4 

(70 


18 

18 

18 

4-5 

(%) 

15 

1 1 

6 

12 

> 5 

(70 

4 

3 

3 

3 

Average At prop ( 1 0 5 y r ) 

2.6 

2.4 


2.4 
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Table 5 

Distributions of SNR Radii 


Rj Bin 
(pc) 


A 

Model 

B 

C 

D 

a) RADII AT DEATH 

25 - 30 (%) 

1 

5 

8 

1 

30 - 35 

(%) 

6 

8 

17 

13 

35 - 40 

(%) 

4 

7 

9 

8 

40 - 45 

(%) 

20 

24 

21 

26 

45 - 50 

(%) 

46 

37 

35 

34 

50 - 55 

(%) 

5 

5 

4 

3 

55 - 60 

(%) 

19 

15 

7 

15 

Average Rj 

at death (pc) 

47 

45 

42 

45 

b) RADII OF ONLY SNRS WHICH MADE INTERSECTIONS, AT t int 


0* 

(%) 

31 

33 

41 

36 

17.6 

(%) 

31 

29 

27 

33 

22.2 

(%) 

7 

8 

7 

6 

27.7 

(%) 

10 

10 

8 

10 

30 - 35 

(%) 

13 

14 

12 

11 

35 - 40 

(%) 

8 

6 

4 

4 

> 40 

(7c) 

0 

0 

0 

0 

Average R { at t jn( (pc) 

17 

16 

14 

15 

c) RADII OF ONLY SNRS WHICH INITIATED SHOCK SYSTEMS, 


AT FIRST t „ 

cav 

i7.6 m 

10 

11 

12 

12 

22.2 

(%) 

11 

7 

8 

7 

27.7 

(%) 

6 

11 

20 

14 

30 - 35 

(%) 

6 

6 

4 

6 

35 - 40 

(%) 

12 

12 

16 

13 

40 - 45 

(%) 

40 

40 

27 

38 

45 - 50 

m 

15 

11 

12 

10 

Average R. at First t cav (pc) 

36 

36 

35 

10 

♦This percentage of SNe occurred within preexisting cavities. 
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Table 6 

Distributions of Effective Volumes of SNRs 


V eff Bin 
{10 s pc 3 } 


A 

Model* 

B 

D 

< 5 

(%) 

... 44 

32 

29 I 

5-10 

(%) 

25 

31 

31 i 

10 - 15 

(%) 

12 

13 

1 6 

15-20 

(%) 

. . . 7 

6 

10 

20 - 25 

(%) 

4 

8 

7 

25 - 30 

(%) 

"t 

... D 

4 

3 

30 - 40 

(%> 

2 

2 

2 

40 - 50 

(%) 

. . . «1 

0 

0 

> 50 

(%) 

2 

3 

3 

Average V cff 

( 10 5 pc 3 ) 

18 

31 

35 

*Data for are not available for Model C. 
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FIGURE CAPTIONS 


FIG. 1. — An Idealized shock system consisting of a long cylindrical 

tunnel and a new SN inside. The inset shows an enlarged 

view of the wall. Gas parameters (see Table 2) are labeled 

for tunnels, walls, and ambient ISM. For a He abundance 

-24 "3 

of 10%, p 0 ~ 1-5 x 10 g cm and P w ~ 2 P c , 

FIG. 2. — The evolution of (a) pressures P^ and (b) radii R_^ for a 

pair of interacting SNRs, as it would be calculated within 

the numerical simulation. One SNR appears at time zero, 

represented by the lower solid curve in (a) and the upper 

5 

curve in (b). The second SNR appears after 4 x 10 yrs at 
a distance of 60 pc from the first explosion, intersects at 
t. nt , and completes a connection at t • The cavity shock 
pressure drop P cav iu the newer cavity, which occurs at , 

is calculated as explained in the text and is equal to the 
vertical distance between the end points of the dashed line. 

At t cav the older cavity receives its pressure impulse, 
defined by the requirement that thereafter the two cavities 
are in pressure equilibrium with each other. After t ccm + 
At^^ the joint cavity expands uniformly to equilibrium 
with P ISM . 

FIG. 3. — The volume fraction in tunnels, averaged over Z, as a 

function of time for four simulated models with average 
porosities q as shown. Each point represents the fraction 
after ten time steps, and the number of SNe in the test space 
during that period is shown in the lower histograms. A burst 
of SNe is directly related to a surge in fraction, and the 
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isolated SNR lifetime and effective lifetime (from Table 3) 
are shown as an aid in interpreting the duration of the 
surges. The lower histograms appear different in all runs 
except C and D, although the input functions were the same, 
because ten time steps is a different fraction of the total 
time in each run. 

FIG. 4. — The quasi-equilibrium tunnel fraction averaged over Z and 

time, as a function of Z-averaged porosity for several runs. 
The dashed line shows the relation expected without SNR 
interaction. 

FIG. 5. — The tunnel fraction f(Z) and porosity q(Z) as functions of 

height above the plane, at the time a typical run (Model B) 

ended. q(Z) is an exponential with scale height 90 pc. 

FIG. 6. — The number of SNRs %<rith diameter less than D as a function 

2.5 

of D. The line labeled D corresponds to adiabatic 

blast waves and is the same as that drawn by Clark and 

Caswell (1976) to fit their sample of young "radio" SNRs. 

The points on the right are (unsealed) results for "X-ray" 

SNRs from the four simulation runs, not from observations. 

The dashed line (broken at the middle) is an extrapolation 

from <-he "X-ray" points down to the "radio" points; the 

3 

former must be scaled by a factor of 10 , roughly the volume 

ratio of galactic spiral arms to test space, in order to join 

the two halves of the line, showing that there is no basic 

inconsistency between radio observations and the simulation 

4 

results. The line labeled D illustrates the steepness of the 
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intermediate-diameter !, X-ray" curves compared to the slope 
expected for SNR shells expanding according to momentum 
conservation „ 
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